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Abstract 

We present a reduced basis approach to solve the convected Helmholtz 
equation with several physical parameters. Physical parameters charac¬ 
terize the aeroacoustic wave propagation in terms of the wave and Mach 
numbers. We compute solutions for various combinations of parameters 
and spend a lot of time to figure out the desired set of parameters. The re¬ 
duced basis method saves the computational effort by using the Galerkin 
projection, a posteriori error estimator, and greedy algorithm. Here, we 
propose an efficient a posteriori error estimator based on the primal norm. 
Numerical experiments demonstrate the good performance and effectivity 
of the proposed error estimator. 


1 Introduction 

Many applications such as estimation of radar cross section [8], heat transfer 
phenomena with high Peclet number [18) . propagation of wave acoustics, and 
so on in physics and engineering, are described by partial differential equations 
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(PDEs) with proper boundary conditions, 


Cu = 

f, 

in 17, 

(la) 

Bu = 

9, 

on 917, 

(lb) 


where C and B are operators for functions on and dVL, respectively, and 17 is 
the domain of the problem with a boundary 917. 

To reflect the physical and geometrical changes and evaluate their effect on 
the result, we introduce input parameters and outputs of interest which are just 
parameters and functional values of solutions. Input parameters are divided 
into physical and domain parameters. The change of physical parameters such 
as density, porosity, frequency, absorption coefficient, flow rate, etc. depending 
on the problem, corresponds to the change of the operators from £, B, /, g to 
£^, /^, 5 ^, where the subscript g denotes the parameter. The deformation 

of the geometrical configuration caused by varying of domain parameters m 
is also studied by the geometric parametric variation [23) of the domain and 
its boundary denoted by 17^ and dftfj,. The dependence on parameters leads us 
into a parametrized partial differential equation (P^DE) from (IT|) 




in 17^, 

(2a) 


9 ( 1 .) 

on 917^, 

(2b) 


and a functional value = l{u^), where I is a functional of interest, and 
denotes the solution depending on the parameter. The output can be statistical 
when the input is stochastic as treated in [113 E]. 

Among many aspects to view the P^DE, there are two main contexts, so 
called the real time and many query contexts [20l [23] , to be considered crucial 
at least in computational engineering. The former is found in parameter estima¬ 
tion or control problem, interpreted as “deployed” or “in the field” or “embed¬ 
ded.” That is, the parameter must be estimated rapidly “on site”. Meanwhile, 
the latter is pursued in design optimization or multi-scale simulation. The state 
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equations should be solved for many parameters [mini US] in the optimization 
problem, and many calculations of small scale problems are required to predict 
the macro scale properties in the multi-scale simulation. Following these con¬ 
texts, the P^DE should be solved rapidly without severe loss of reliability, that 
is, while keeping the almost same order of approximation, the evaluation must 
be done as soon as possible. 

According to [55] [53] , we can regard the set of solutions generated by pa¬ 
rameters in a parameter domain as a smooth low-dimensional manifold in the 
approximate space. The reduced basis method (RBM) is based on the low order 
approximation of the manifold owing to the low dimensionality of the solution 
manifold. Under some sufficient assumptions, the computational task of the 
P^DE is decomposed into the off-line and on-line stages. In the parameter in¬ 
dependent off-line stage, a heavy computation is done to generate a reduced 
basis. In the on-line stage, the computation for new parameter is performed 
by the Galerkin projection into the reduced basis space. The marginal number 
of computations gets important since it says about the minimum number of 
computations by the usual method to exceed the total cost of the RBM due to 
the off-line stage. Because of the invention of a posteriori error estimators, rig¬ 
orous error bounds [7] for outputs of interest, and effective sampling strategies, 
the RBM evaluates the reliable output for many combinations of parameters 
in high dimensional parameter space rapidly, which means that the marginal 
number of computation gets smaller. The reliability of the result by the RBM 
is guaranteed by theoretical results in mmm- 

One sufficient assumption to ensure the decomposition of the computational 
task is the affine dependence [22123] of the P^DE, i.e., the related forms are 
expressed by the linear combination of parameter independent forms with pa¬ 
rameter dependent coefficients. Under this assumption, the error bound has 
many terms depending only on the dimension of the approximate space which 
are independent of the parameters, and computed during the off-line stage. This 
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is good point, but there are two bottlenecks in the computational point of view. 
Firstly, the error bound formula is very sensitive to round-off errors, which may 
show a little bit large discrepancy between the a posteriori error bound and 
the on-line efficient formula. Secondly, the RBM is intrusive, which means that 
computation of the solutions requires intervening the matrices assembly routines 
of the code. To remove the intervention, one can use the empirical interpolation 
method [mnnn] which separates the parameter and the space variable of the 
affine coefficients. 

In this paper, we describe the propagation of acoustic waves in a subsonic 
uniform flow by the time harmonic linearized Euler equation and transform it to 
a convected Helmholtz equation for the pressure field in Section [5] The problem 
of the convected Helmholtz equation is well posed when appropriate boundary 
conditions are imposed, see [5] for details. We present a RBM for solving the 
convected Helmholtz equation with these two physical parameters. Physical 
parameters are the Mach and wave numbers, which are the ratio of the mean flow 
velocity and frequency to the sonic speed in the flow, respectively. The outline 
of the RBM is presented with the greedy algorithm in the pseudo code style in 
Section [3l We present numerical simulations by varying physical parameters in 
Section m Finally, several conclusions and future works are addressed on the 
convected Helmholtz equation with several parameters. 

2 Convected Helmholtz Equations 

2.1 Bounded domain 

We consider compressible flows induced in a uniform subsonic flow in the di¬ 
rection of xi with Mach number 0 < M < 1 for (xi,X 2 ) G Assume that 
perturbations in the density p, the pressure p and all components of the velocity 
vector u := (ui, U 2 ) are small, and all sources and initial disturbances bounded 
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to the rectangular domain 


ft = [—oi,ai] X [— 02 , 02 ], 01,02 > 0. 

After nondimensionalizing appropriately, the flow is governed by the linearized 
Euler equation 


Dtu -f Vp 

= 0, 

(3a) 

Dtp -1- V • u 

= 0, 

(3b) 

Dtp -1- V • u 

= 0, 

(3c) 


where Dt = dt + Mdx^ is the convected derivative or the material derivative 
in the direction of (M, 0), see [H] for a detailed derivation of the equation. 
Applying Dt to (l3bl) and V- to (I3a1) . and subtracting between them yields the 
convected wave equation 

idt+Mdxtfp = Ap. (4) 


The Fourier transform of (|4]) with respect to time t gives the convected Helmholtz 
equation 

Usually, we impose a proper boundary condition to solve ([S]). For notational 
convenience, p is used instead of p, then after enforcing a general function / on 
the right-hand side of (O, it takes the following divergence form 


V\MVp + Bp) + k^p = f, 


( 6 ) 



1 

0 

(M 

1 

1—I 

1_ 


2ikM 

M = 

0 I 

, B = 

0 


In one parameter problem of (ED, the wave number k changes under a fixed 
Mach number M. Both M and k varies in their domains of parameters in two 
parameters problem. In this paper, we consider k or (fc, M) as a parameter 
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/i. The variational problem of the convected Helmholtz equation ([S]) is to find 
p G such that for given 0 < M < 1 and fc > 0, 

— j + Bp) ■ Vv dx + j pvdx= f fvdx, for all z; G (7) 

where n is the outer normal vector. 


2.2 Unbounded domain 

As in m, we use the following notations 

rib = {{xi.x^) G : a;_ < xi < x+, —d < X 2 < d}\B, 

^PML _ ^(^xi,X2) : X- — L < xi < X-, —d < X2 < d} 

^PML _ ^(^xi,X 2 ) GM)^ : x+< xi < x++ L,—d < X 2 < d}, 

where B is an obstacle such as a circular or elliptical hole, G M and L > 0. 
From [IH], a PML formulation for the convected Helmholtz equation is 


(l-JlP) o(*i) — + 


d 


ikM 


dxi 1 - M2 


^ d'^p 


( 8 ) 


(9) 


in fi = Hf, U U Here, the damping function is of the form 

-iuj 

o:{xi) = ——^ 

— lU! + ^(Xl) 

a{xi) = <Jo{{xi - X-)'^X{x_-L,x_){xi) + {xi - X+)'^X{x+,x++L){xi))llO) 

where uo is a parameter for the magnitude of damping and xa(x) is the charac¬ 
teristic function on the set A C M. See [12] for other type of the PML condition. 
The divergence form of ([8]) is 


dcy. 

V-(AfaVp-f Hap)-b ( - ifcM^— )p =/, in H, 


( 11 ) 


where 


Ma = 


(l-M2)a 0 


a 


-1 


B^ = 


2ikMa 

0 


2 fc2(a ^ — aM2) 


k^ = 


1 - M2 


6 










Note that m is the same as ([HI) in the region ili, from the definition of the 
damping function. And the variational form of dill) is 

— / {M-a^p-\-BaP)-Vv dx — ikM [ — — pvdx+k'^ [ pvdx= [ fvdx, (12) 

Jh Jh dxi Jq Jq 

for all V S see [H] for details. The equation (fT^ is also the same as (O 

when the support of the test function v is in Vlh. 

3 Reduced Basis Method 

In general, the RBM constructs the reduced basis using the greedy algorithm 
and precompute the parameter independent parts of matrices at the off-line 
stage. We assemble the matrices using the coefficients at new parameter, solve 
the system and compute the output at the on-line stage. In the whole process, 
we restrict the approximate space to the much smaller subspace chosen by the 
greedy algorithm and discard the unnecessary modes during the calculation of 
the basis. The a posteriori estimator measures errors of approximation and is 
the key to the model order reduction [HElEllTlElIIllEllEslEl]. The former 
is given by the property of the approximate space and chosen under the proper 
assumption. The latter depends on the reduced basis subspace. 

3.1 Primal and Dual Problems 

Let X be Hq{V,) with the inner product (-jOjc its associated norm IHIjf. 
Let ^ be a parameter selected from a certain parameter set V. We solve the 
parametrized variational form for (I5|) such as 

a {p{fj,), v; p) = f (v; p), for all £ A 

where and f{-',p) are bilinear and linear forms depending on the pa¬ 

rameter vector /r, respectively. We evaluate the quantity of interest s{p) as the 
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value of a linear functional I G X' at the solution p{fj,) 


s{n) = l{p{p)]p). 

The finite dimensional approximation p^{p) of p(/r) in a smaller function space 
X-^ C X oi dimension Af satishes 

a {p'^ip), v; p) = f{v, p), for all v G X^, (13) 

and its quantity of interest s^{p) is 

s^{p) = ^{p^{p)-.p)- 

The approximate solution p^{p) of (1131) is the truth approximation, which is 
accurate enough for all parameters p G T). To claim the accuracy, we must 
choose a very large Af and thus need to solve a large sparse matrix system of 
algebraic equations. 

In the RBM, we want to make a much smaller space than the approxi¬ 
mate space X^. The space Xjv is called a reduced basis, spanned by the linearly 
independent approximate solutions {p^{Pj)}^^^, he., Xm = Span({p^(/rj)}^^). 
For the user-chosen parameter p G V , the reduced basis approximation pAr(/i) G 
X]Y is obtained by the Galerkin projection, 

a{PNip),v;p)=f{v;p), for all u G X^, (14) 

and its quantity of interest SAr(/r) is 

sn{p) = l{pNip);p)- 

Note that the reduced basis space X^r of dimension N is much smaller than the 
finite approximate space X^ of dimension Af. 

To improve the order of convergence of output, i.e., quantity of interest s{p), 
we introduce the dual problem of the primal problem (1131) : find {p) G X^ 
such that 

a [v,w^(p)] p) =—l{v; p), for all u G X^. (15) 
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Its reduced basis approximation of is also defined by the Galerkin 

projection 

a {v,WN{fi); ^) = —l(v; fi), for all G (16) 

Formally, the error and residual relations of the primal problem are written as 


eP(/x) = - pNifi), 

r^(v\p) = f{v]p) - a{pN{p),v\p) = a{p^{p),v]p) - a{pN{p),v,p) 

= a {e^(fi), V, p) , for all ri € Xjv, 


and those of the dual problem are expressed by 
e'i(^) = w^{p)-WNip), 

r‘^{v\p) = -l{v;p)-a{v,WN{p)-,p)=a{v,w'^{p)-,p)-a{v,WN{p);p) 

= a (w, e'^(/i);/i), for all u S Xjv. 

We call eP(/r), e'^{p) and r‘^(-;/x) the primal error, the primal residual, 

the dual error and the dual residual, respectively. As in m Section 2], the dual 
corrected output (p) is defined by 

s^j^ip) = l{pN{p)\p) - r^{wN{p)-,p)- (17) 

Then the error {p) — s^{p) is expressed in terms of the dual residual of the 
primal error, 


s^{p) - s^^ip) 


+r^{wN{p)\p) 

l{e^{p)]p) + a {e^{p),WN{p);p) 

-a (eP(^), w^{p)\p) + a {e^{p),WN{p)\p) 
-a {e^{p),e‘^{p);p) 

-rd(eP(/r);/r), 


and it is bounded by the norms of the primal error and the dual residual 


s^{p) - s'^^ip) < ||r‘^(-;Ai)||^,||eP(/r)||^, 
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where the dual norm ||/||j 5 ^/ of any linear functional I G X' is defined in the usual 
sense: 


lx- 


= sup 


l{v) 


vex ||v|lx 

Note that there is improvement in the convergence by the solution of a dual 
problem, see [531 Section 11] for more details. To treat the non-coercive problem, 
we may assume that the bilinear form of the system satisfies an inf-sup condition. 
The non-zero inf-sup stability constant /3(/r) of a (•, •; 


Iblljf <sup for all p ex, 

/3(p) Ibllx 


makes it possible to bound the norm of the primal error by the dual norm of 
the primal residual 


l|eP(M)llx < 


1 |a(eP(p),-c;p)| 

/3(m) .ex Ibllx 


1 \r^{v;f^)\ 

/3(M)®ex Ibllx 






We can also bound the primal output error by the norms of the output and the 
primal residual. 


s^(p)-s^(/r)| = beP(p);p)|<||l(-;/r)||x-||eP(/r)||^ 

- ^ll^(•;A^)llx'lb"(•;A^)llxo 


and the dual corrected output error by the norms of the primal and dual resid¬ 
uals. 




< 


1 


|r'i(.;/r)|L,||rP(.;/r)|U,. 


The Riesz representation eP(/r.) G X of the primal residual rP(-;/r) G X' such 


that 


{e^(fi), v)^ = rP {v, fi), for all u G X, (18) 

is very useful for the computation of the error estimators of the off-line and 
on-line stages in the RBM. From the definition of the dual norm, we obtain 
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that 


lk‘'(-;/^)llx' = sup 

x^X 


rP{v,fi) 

llullx 


sup 


Hull;, 


\\en^^)\\x■ 


3.2 Matrices and computational costs 

Let {Cj}^i be the orthonomalized reduced basis for Xx = Span({p^(^j)}^^), 
where are parameters selected from some sampling strategy. Then we 

may expand the solution pAr(/i) for the parameter /i S in terms of this reduced 
basis 

N 

PN{f^) = ( 19 ) 

i=i 

Inserting this expansion into m and applying a test function gives us the 
fc-th row of the TV-dimensional system of equations: for k = 1,..., 

N 


E^T(M)a(Ci,Cfc;M) =/(Cfc;M), (20) 

J=1 

with the following output 


N 

sn{p) = 

i=i 


( 21 ) 


Denote by $ the matrix consisting of the reduced basis as its column 


vectors, 


^’ = [Cl • • • Cn] 


then due to the orthonormal property of the reduced basis in X, it satisfies 


($*,$);, = l;v, 

where $* is the Hermitian of $, and Iat is the identity matrix of order N. Note 
that the inner product of X is extended to the matrices of order N. Using the 
coefficient vector C(/i) whose components are C; (m)i we may rewrite ([19]), (l20|) 
and m as follows: 

Pn{p) = $C(m), ^M(/x)$C(m) = sxip) = L{n)*<^^{p), 
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where A(/i), F{^) and i(/i) are the matrices representing /(s/^) and 

Here, the matrices satisfy the following properties: for any p = 

V = ^rj G X]\j, 

{<i>r])*A{p){<i>^)=a{p,v;p), {^^)*F{p) = f{p; p), L{p)*^^^l{p-p).{22) 

We can define the Riesz representation e^{p) of the primal residual in (I18p 
explicitly, 

ePip)^Z-^[F{p)-Aip)<^>a^^)], 
where Z is the matrix due to the inner product such that 

($? 7 )*Z($ 5 ) = (p, w)^, for all p = $^, u = $77 e Xn, 

and Z* = Z from the property of the inner product in X. Then the square 
norm ||eP(p)||^ is 

\\e^p)\\l = F(p)*Z-iF(p)-F(p)*Z-iH(p)$e(p) 

-m^^)rA{prz-^F{p) + mf^)rA{prz~^A{p)<^af^). 

Similar to the primal problem, let {C 7 be the reduced basis from solutions 
{p) of (fTKll for the same parameters {pj}^^i- Let $d and ^d{p) be the matrix 
for the reduced basis and the coefficient vector for wn{p) of (fTBl) . then 

we can write 

wn{p) = $d^d(p), $d-4(p)$d&(p) = -^lL{p). 

Using (1^ . the dual corrected output dm becomes 

s^j^{p) = - {^dUp)nF{p) - H(p)$e(p)). 

And let e^{p) be the Riesz representation of the dual residual 
(u, e‘^(p))^ = p), for all v G X. 

It is expressed as 

e\p) = Z-^[-L{p) - AipY^diM. 
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with the square norm ||e‘^(/i)||^ is 


In the RBM, it is very crucial to assume that all the related forms may 
be expressed as the linear combinations of parameter independent forms with 

parameter dependent coefficients, or they may be affine in the parameter: 

/ 

a {p, V; p) = Em=l ^a,m{p)am (p, v), 

< /(w;m) = Em=l 0/,m(A^)/m(^)i (^3) 

KP\ P) = EmLl ^l,m{p)lm{p)- 

Here, am {■,■), fm{') and lm(,-) are parameter independent forms. Clearly, 
Qa,m{p), Qf,m{p) and 0z,m(/i) are parameter dependent coefficients. This as¬ 
sumption enables us to realize an efficient off-line and on-line splitting during 
the computational procedure. The above is expressed in matrices 

Ma Mf Ml 

A{p) = ^ ^ ^ip) ~ ^ ^ ^f,7nip)^nii ^ip) ~ ^ ^ ^l,Tnip)^7ni 

m—1 m—1 m—1 

where Am, Fm and Lm are the matrices representing am ('jOi /m(’) and Imi')- 
When the related forms are affine as in (l23ll . the approximate system (l20ll be¬ 
comes 

M„ N Mf 

EE Qa,mip)^jip)0‘m iCj , Ck) — ^ ^ ^ f,mip) fmiCk) ■ 

m—1 j—1 m—1 

Let fm and dmj be the Riesz representations of /m(’) and am iCji ') such that 

(/m,w)^ =/m(u), idm,j,v)j^ = amiCj,v), for all w £ X. (24) 

Then we have the following representation of ePip), 

Mf Ma. N 

s^ip) — 'y ^ 0/,m(li)/m ~ 'y ^ ^^^a,mip)^jip)dm,j, 

m—1 m—1 j—1 
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and its norm may be expressed as 

Mf 

m,n—l 

Mf Ma N 

— Q/,m(/^)Qa,K(/i)$j(/^) f/mi ^n,j) __ 

m—1 n—1 j—1 
Mf Ma AT 

-EEE 

n=l m=l j = l 
Ma AT 

+ E E 0a,m(/^)0a,n (/^)Cj (/^)Cfc (/^) (^m,j j ^n,j)x’ 

m,n—l j,k—l 

which is independent of A/” after off-line computations of the Af dependent quan¬ 
tities /m and a^j with the inner products (/m,/«)x> (/m)a«j)x> {^m,j,fn)x 
and (am,j, ^n,j)x- The number of operations to evaluate ||eP(/r)||j^, or the com¬ 
putational cost C(||eP(/x)||jj.) is 

C(||eP(Ai)||^) = 3Mf + 8MfMaN + 5M2N^ 

where the operational costs of addition, subtraction, multiplication and square 
root are assumed to be of the same order. The coefficients ^j(/r) of PNifJ') are 
obtained after solving the reduced system (m of dimension N, whose cost is 
denoted by Cx- In many cases, Cjv rnay not be of order due to the lack of 
the sparsity of the reduced system (HI- 

At the off-line stage, we solve approximate solutions satisfying (TT^ to form 
the reduced basis and orthonomalize them. Using the reduced basis we 

need to compute the Riesz representations ot Mf + MaN forms and the inner 
products of -|- 2MfMaN + pairs. Thus the computational cost Cog 

at the off-line stage is 

Coe = NCx + Cq + {Mf + MaN)CR + [M] + 2MfMaN + M^N^) {2N - 1), 

where Cj\f, Cq, and Cr are the computational costs to solve the system (fT^ of 
dimension Af, orthonomalize Xx including a posteriori error estimators, and 
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compute the Riesz representation in respectively. When the system m 
is sparse, is of order A/”^. If the Riesz representation is bounded in then 
Cr is of order N. 


3.3 Error estimator and greedy algorithm 


Examining the bounding formula of the primal error, we may define the following 
error estimator and its effectivity: 


A(/r) 


1 






A(m) 

\\p^{f^) - Pn{Mx' 


where the effectivity quantifies the performance of the error estimator A(^) 
for the reduced basis solution. The stability constant is assumed to be constant 
during the calculations, which causes slight loss of effectivity but still works. 
We judge the current reduced basis approximation is sufficiently accurate if all 
values of the selected error estimator are smaller than the given tolerance. 

In Algorithm [T] (Greedy Algorithm) [6l [7l [9l [22 [24], it starts from the 
selection of the training sample set Ptrain from the parameter space T), the tol¬ 
erance e of the error estimator and the maximum dimension Amax of the reduced 
basis at lines [2111 For the randomly chosen parameter fj. in Strain, we compute 
the solution of (fl^ and normalize it with respect to the inner product (•, Ox of 
X at lines [2Q3 Then we search the next parameter maximizing the error esti¬ 
mator among parameters in Ptrain at line 1 1 21 If the error estimator for the new 
parameter is smaller than the tolerance or the dimension of the basis system is 
over the maximum dimension, then the process stops at line \14\ Otherwise, we 
compute the solution, orthonormalize the new basis including the previous ones, 
construct the error estimator, find the next parameter maximizing the error es¬ 
timator, evaluate the error estimator for the new parameter, and examine the 
result to decide whether to stop the process at lines [T3H7PI and line [HI Using 
the final reduced basis, we can compute the solution and the residuum for new 
parameter at lines [2411271 Algorithm [T] shows these procedure for the problem 
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m in the pseudo code style. 
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Algorithm 1 Greedy Algorithm 
1: procedure Initialization: 

2: Construct the training sample set Strain 

3: Specify a tolerance e as stopping criteria 

4: Choose the maximal dimension of the reduced basis space 

5: end procedure 

6: procedure Offline procedure: 

7: Choose the first parameters fXi randomly 

8: Compute the snapshot p(/ri) of ([S]) 

9: Set Xi = {p(Ati)} 

10: Orthonormalize Xi 

11: Construct the residuum Ai 

12: fi 2 = avg max Ai(/r) 

13: i = 2 

14: while Ai_i(/ii) > e and i < Amax do 

15: Compute the snapshot of ([5]) 

16 : Set Xi = Ai_i U {p(/ri)} 

17: Orthonormalize Xi 

18: Construct the residuum A^ 

19: /Xi +1 = arg max Ai{fi) 

20: i = i + l 

21: end while 

22: end procedure 

23: procedure Online procedure: 

24: Choose new parameter fj, 

25: Compute the coefficients of the reduced system 

26: Determine the solution using the reduced basis 

27: Compute the residuum of the solution 

28: end procedure 
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4 Numerical Results 


The fundamental solution of the convected Helmholtz equation generated from 
the point source at the origin is 


^(xi,X2) = 


i ( 1 ) / kyjxl + {I - M‘^)xl 


1 -M2 


exp I —i 


kMxi 


1-M2 


where (z) is a Hankel function. 

Let K he a, triangular element consisting of three vertices {xi,yi), {x 2 ,y 2 ) 
and (xsji/a). Obviously, it belongs to the triangulation T of H and \K\ denotes 
its area. The affine transformation from the reference triangle K to the 
triangle K is defined by 


:x^ T^{x) = B^x + b^, 


where 


II 

X 2 — Xi X 3 — Xi 

II 

Xi 


_ y2 -yi 2/3 - yi _ 


_ yi _ 


We use the PI conforming finite element basis function. Then from ([7]), the 
local system at the element K satisfies, for a local solution vector , 

+ M^ +C^)p^ = f^, 

where the stiffness , mass and convection matrices are 

^^af{M)Su = =-ikM - B^^C2), 

i=i 

where the parameter independent parts are 



1 


-1 

0 


2 

-1 

-1 


1 

0 

-1 


Si = 

-1 


1 

0 

, 52 = 

-1 

0 

1 

, 53 = 

0 

0 

0 



0 


0 

0 


-1 

1 

0 


-1 

0 

1 



2 

1 

1 



-1 

-1 

-1 


-1 

-1 

- 

1 

M = 

1 

2 

1 


Cl = 

1 

1 

1 

, C 2 = 

0 

0 

0 


1 

1 

2 



0 

0 

0 


1 

1 

1 
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with parameter dependent coefficients af {I = 1, 2, 3) 


if Z = 1 


af (M) = X - ^ 21^22 + if ; = 2, 


4\K\ 


{B^,)\b^,Y -M\B^,y 


a I = 3. 


Similarly, we can express the force vector after imposing appropriate bound¬ 
ary conditions for the boundary integrals in © such that simple Dirichlet con¬ 
dition. We assemble these local systems into the global system and solve the 
problem for the given parameter. For the PML case, we can also derive a similar 
afhne system of the linear combination of parameter independent matrices and 
parameter dependent coefficients. 

The computational cost for one realization of uncertain parameters in the 
problem by the Galerkin method is lower than the total computational cost 
including the off-line cost CoS and the on-line cost Con by the RBM, but if we 
want to solve the problem with many different realizations of parameters, the 
reduced basis allows us to reduce the total computational cost. Let n be the 
number of computations due to realization of parameters. Then the RBM is 
profitable when uCg > CoS + nCon, where Cg is the computational cost of the 
Galerkin method. In short, the profit by the RBM occurs whenever 


n > 


Coff 


Cg — Co 


(25) 


holds. We call the smallest integer n* satisfying the inequality (Ell) as the 
marginal number of the RBM, which indicates the minimum number of compu¬ 
tations to get the benefit of the RBM in the aspect of the total computational 
cost. Usually, the computational costs are measured in seconds. 
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Figure 1: Bounded domain excluding a circular hole. 

4.1 Bounded Domain 

The bounded domain is a box [—1,1] x [—1,1] except a circular hole of radius 0.3 
and center at the origin as shown in Figure [TJ We choose between h = 0.03 
and h = 0.025027 as the maximum diameter of elements in the mesh, which 
is called by the mesh size. The first choice generates 4800 vertices and 9250 
elements, while the latter does 10841 vertices and 21157 elements in the domain 
by Gmsh [TU] . 

4.1.1 One Parameter of Wavenumber 

We use the training sample set Ptrain = consisting of an Ni terms of an 

arithmetic progression sequence from to fcmax, where Ni, fcmin and fcmax 
are the number of samples, lower and upper bounds of fc in ((S]), respectively. We 
take 40 samples (iVi = 40) from the interval between k^in = 2 and fc^ax = 5. 
Numerical computations are done for M = 0.3 and M = 0.4 in the mesh of 
h = 0.03. 

For M = 0.3, Figure [2] comes from the residuum columns of Table [T] and 
illustrates the evolution of the residuum as the dimension N of the reduced basis 
space increases. It shows very fast decrease of residuum after the dimension 
exceeds 16. We report that the computational costs at the on-line stage are 
Con = 0.0348 for M = 0.3 and Con = 0.0408 for M = 0.4. Compared to 
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the computational costs of one computation by the Galerkin method, those at 
the on-line stage are 180 and 150 times shorter for M = 0.3 and M = 0.4, 
respectively, which can be calculated by taking ratios of Cq in Table [T] to Con- 
The marginal number n* increases as the dimension of the reduced basis space 
does. For instance, the marginal numbers n* are 876 and 871 for M = 0.3 and 
0.4, respectively, at the 28 dimensional reduced basis space. 



Figure 2: Convergence in function of the dimension of the reduced basis space for a 
parameter k under M — 0.3. 

4.1.2 Two Parameters of Wave and Mach Numbers 

Let Utrain = ^ made of the product of an Ni terms of an 

arithmetic progression sequence from fcmin to fcmax, and an N 2 terms of an 
arithmetic progression sequence from Mmin to Mmax- Here A^i, N 2 , /cmin, fcmax, 
Mmin and Mmax are numbers of samples in wave number k and Mach number 
M in dSD, lower bounds and upper bounds of them, respectively. We set = 
N 2 = 10, /cmax = 12, fcmin = 8, M^ax = 0.4 and Mynin = 0.2 in the mesh of size 
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Convergence 



Figure 3: Convergence in function of the dimension of the reduced basis space for a 
parameter k under M — 0.4. 


Table 1: Residuum, computational costs of the off-line stage and one Galerkin solution 
for M = 0.3 and M = 0.4 according to the dimension of the reduced basis space. 



M = 0.3 

M = 0.4 

Dimension 

Residuum 

Coff 

Cg 

Residuum 

CoS 

Cg 

4 

2.7627 

789.4 

6.37705 

2.8557 

793.9 

6.42678 

8 

2.3673 X 10-2 

1581.1 

6.36272 

1.0615 X 10-2 

1590.0 

6.41576 

12 

1.3459 X 10-3 

2370.7 

6.37737 

3.1121 X 10-3 

2384.2 

6.37039 

16 

4.9738 X 10-3 

3163.7 

6.37565 

2.4937 X 10-'‘ 

3181.8 

6.38359 

20 

2.4432 X 10-3 

3958.1 

6.38165 

4.4107 X 10-3 

3979.2 

6.42344 

24 

2.9129 X 10-“ 

4750.5 

6.38063 

2.6859 X 10-“ 

4778.2 

6.43455 

28 

8.0662 X 10-“ 

5542.6 

6.37059 

5.4560 X 10-“ 

5575.4 

6.44458 
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h = 0.025. 

The computational costs at offline CoS, online Con and one full Galerkin 
method Cq are 23281, 0.1099 and 25.538, respectivley. We see that the com¬ 
putational benefit of the RBM occurs when the computations are more than or 
equal to the marginal number n* = 916. We also see that the computational 
cost at the on-line stage is 230 times shorter than that of one computation by 
the Galerkin method, where the number of the speed up comes from the ratio 
of the Galerkin cost Cq to the on-line cost. This is very promising aspect of 
the RBM such that the speed up makes it possible to apply the RBM to the 
practical problems under many and fast computational loads. 

Figure 0] shows the real part of the solution by the reduced basis of dimen¬ 
sion = 10 and the absolute error between the RBM solution and the exact 
solution for fixed parameters M = 0.3 and k = 10. The errors between the 
RBM solution and the exact one are 0.0278 in 0.0223 in and 

0.0320 in 

4.2 Unbounded Domain 

The duct in Figure [B] has an elliptical hole whose major and minor axes are 
a = 0.3 and b = 0.25, and center is at the origin. We set X- = —1, x+ = 1, 
L = 1 and tJo = 15 for the damping function a{x) in We generate meshes 
for Cl of mesh size h = 0.0381 by Gmsh, which has 16907 nodes and 33262 
elements. We treat the wave and Mach numbers as parameters. We use 16 
training samples among [8,12] x [0.2, 0.4] and choose 10 basis from them. The 
computational costs at offline Coir, online Con and one full Galerkin method Cq 
are 4899, 0.15763 and 188.2764, respectivley. The marginal number is n* = 27 
and the computational speed by the on-line stage is at least 1,100 faster than 
that by the usual Galerkin method. 

The errors between the 10 dimensional RBM solution and the exact one are 
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Figure 4: Real part of the numerical solution by the 10 dimensional RBM chosen 
from 100 samples among [8,12] x [0.2, 0.4] when M = 0.3 and k = 10. 
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Absolute Error 
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Figure 5: Absolute error of the numerical solution by the 10 dimensional RBM chosen 
from 100 samples among [8,12] x [0.2, 0.4] when M = 0.3 and k = 10. 




Figure 6: Bounded domain fit and PML domain U 
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0.0682 in 0.0248 in L^{n), and 0.4012 in H^{n). The H^{n) error is 

higher than that for the bounded domain, which is caused by the small number 
of training samples. 

^ Approximate Solution 

0.8 

0.6 

0.4 

0.2 
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- 0.2 

- 0.4 

- 0.6 

- 0.8 

-1 

-1 - 0.5 0 0.5 1 

Figure 7: Real part of the numerical solution by the 10 dimensional RBM selected 
from 16 samples among [8,12] x [0.2, 0.4] in fib when M = 0.3 and k = 10. 



5 Conclusion 

We test the RBM for the convected Helmholtz equation. The physical parame¬ 
ters are expressed as coefficients of the equation. After these tests, we confirm 
that the RBM works well and gives us the benefit of fast computation at least 
100 times than the usual computational method does. In the implementation, 
we use the error estimator based on the primal norm of the error. 
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Absolute Error 



Figure 8: Absolute error of the numerical solution by the 10 dimensional RBM selected 
from 16 samples among [8,12] x [0.2, 0.4] in Qi, when M — 0.3 and k = 10. 

References 

[1] M. Barrault, Y. Maday, N. C. Nguyen, and A. T. Patera, An ‘empirical 
interpolation’ method: application to ejjicient reduced-basis discretization 
of partial differential equations, C. R. Acad. Sci. Paris Ser. I 339(2004), 
no.9, 667-672. 

[2] E. Becache, A.-S. Bonnet-Ben Dhia, and G. Legendre, Perfectly matched 
layers for the connected Helmholtz equation, SIAM J. Numer. Anal. 
42(2004), no.l, 409-433. 

[3] P. Binev, A. Cohen, W. Dahmen, R. A. DeVore, G. Petrova, and Prz. Woj- 
taszczyk. Convergence rates for greedy algorithms in reduced basis methods, 
SIAM J. Math. Anal. 43(2011), no.3, 1457-1472. 


27 





[4] S. Boyaval, C. Le Bris, T. Lelievre, Y. Maday, N. C. Nguyen, and A. T. 
Patera, Reduced basis techniques for stochastic problems, Arch. Comput. 
Methods Eng. 17(2010), 435-454. 

[5] S. Boyaval, C. Le Bris, Y. Maday, N. C. Nguyen, and A. T. Patera,A 
reduced basis approach for variational problems with stochastic parameters: 
Application to heat conduction with variable Robin coefficient, Comput. 
Methods Appl. Mech. Engrg. 198(2009), 3187-3206. 

[6] A. Buffa, Y. Maday, A. T. Patera, Chr. Prud’Homme, and G. Turinici, 
A priori convergence of the greedy algorithm for the parametrized reduced 
basis method, ESAIM Math. Model. Numer. Anal. 46(2012), no.3, 595-603. 

[7] Y. Chen, J. S. Hesthaven, Y. Maday, and J. Rodriguez, Certified reduced 
basis methods and output bounds for the harmonic Maxwell’s equations, 
SIAM J. Sci. Comput. 32(2010), no.2, 970-996. 

[8] Y. Chen, J. S. Hesthaven, Y. Maday, J. Rodriguez, and X. Zhu, Certified 
reduced basis method for electromagnetic scattering and radar cross section 
estimation, Comput. Methods Appl. Mech. Engrg. 233-236(2012), 92-108, 
2012 . 

[9] M. Drohmann, B. Haasdonk, and M. Ohlberger, Reduced basis approxi¬ 
mation for nonlinear parametrized evolution equations based on empirical 
operator interpolation, SIAM J. Sci. Comput. 34(2012), no.2, A973-A969. 

[10] C. Geuzaine and J. F. Remade, Cmsh: a three-dimensional finite element 
mesh generator with built-in pre- and post-processing facilities, Int. J. Nu¬ 
mer. Meth. Eng. 79(2009), no.11, 1309-1331. 

[11] B. Haasdonk, K. Urban, and B. Weiland, Reduced basis methods 
for parametrized partial differential with stochastic infulences using the 


28 



Karhunen-Loeve expansion, SIAM/ASA J. Uncertain. Quantif. 1(2013), 
79-105. 

[12] F. Q. Hu, On absorbing boundary conditions for linearized Euler equations 
by a perfectly matched layer, J. Comput. Phys. 129(1996), 201-219. 

[13] Chr. Jaggli, L. lapichino, and G. Rozza, An improvement on geometri¬ 
cal parameterizations by transfinite maps, C. R. Acad. Sci. Paris Ser. I, 
352(2014), 263-268. 

[14] Y. Maday, A. T. Patera, and G. Turinici, A priori convergence theory for 
reduced-basis approximations of single-parametric elliptic partial differential 
equations, J. Sci. Comput. 17(2002), 437-446. 

[15] A. Manzoni, A. Cohen, and G. Rozza, Shape optimization for viscous flows 
by reduced basis methods and free-form deformation, Int. J. Numer. Meth. 
Fluids, 70(2012), 646-670. 

[16] F. Negri, A. Manzoni, and G. Rozza, Reduced basis approximation of 
parametrized optimal flow control problems for the stokes equations, Com¬ 
put. Math. Appl. 69(2015), 319-336. 

[17] F. Negri, G. Rozza, A. Manzoni, and A. Quarteroni, Reduced basis method 
for parametrized elliptic optimal control problems, SIAM J. Sci. Comput. 
35(2013), no.5, A2316-A2340. 

[18] P. Pacciarini and G. Rozza, Stabilized reduced basis method for 
parametrized advection-diffusion PDFs, Comput. Methods Appl. Mech. En- 
grg. 274(2014), 1-18, 2014. 

[19] S. Park and I. Sim, Analysis of PML method for stochastic convected 
Helmholtz equation, submitted. 


29 



[20] A. T. Patera and G. Rozza, Reduced basis approximation and a posteriori 
error estimation for parametrized partial differential equations, Version 1.0, 
MIT Pappalardo Graduate Monographs in Mechanical Engineering, 2007. 

[21] J. Pomplun and F. Schmidt,AcceZeroted a posteriori error estimation for 
the redueed basis method with applieation to 3D eleetromagnetic seattering 
problems, SIAM J. Sci. Comput. 32(2010), no.2, 498-520. 

[22] A. Quarteroni, G. Rozza, and A. Manzoni. Certified reduced basis approx¬ 
imation for parametrized partial differential equations and applications, J. 
Math. Ind. I(20II), Art. 3. 

[23] G. Rozza, D. B. P. Huynh, and A. T. Patera, Reduced basis approximation 
and a posteriori error estimation for affinely parametrized elliptic coercive 
partial differential equations. Arch. Comput. Methods Eng. 15(2008), no.3, 
229-275. 

[24] D. Wirtz, D. C. Sorensen, and B. Haasdonk, A posteriori error estimation 
for DEIM reduced nonlinear dynamical systems, SIAM J. Sci. Comput. 
36(2014), no.2, A311-A338. 


30 



